#pragma once

#include "Integrator.h"

/**
 * @file GhostFiller.h
 * @author xingyifan
 * @date 2022-08-31 22:39
 *
 * @description: ghost cell 填充器
 */

template <int Dim>
class GhostFiller
{
public:
    // 默认使用周期边界
    GhostFiller()
    {
        for (int d = 0; d < Dim; d++)
        {
            type_[2 * d] = 2;
            type_[2 * d + 1] = 2;
        }
    }
    GhostFiller(FuncXt g, int *type) : g_(g)
    {
        for (int d = 0; d < Dim; d++)
        {
            type_[2 * d] = type[2 * d];
            type_[2 * d + 1] = type[2 * d + 1];
        }
    }

    /**
     * @file GhostFiller.h
     * @author xingyifan
     * @date 2022-08-31 22:39
     *
     * @description: 使用此函数获得十字形 ghost cell
     */

    Vector operator()(Tensor<Dim> &v)
    {
        // 创建容器，不需要初始化
        Vector ghost(5 * Dim);

        // 循环每一个方向，填充 ghost cell
        for (int d = 0; d < Dim; d++)
        {
            // 获得指定方向上的 ghost cell
            operator()(v, ghost.ptr() + 5 * d, d);
        }

        return ghost;
    }

    /**
     * @file GhostFiller.h
     * @author xingyifan
     * @date 2022-08-31 12:55
     *
     * @description: 在 d 维上获得 ghost cell
     */

    void operator()(Tensor<Dim> &v, Real *ghost, int d)
    {
        // 网格规模
        int n = v.width();

        // 在此方向上移动一个点的位置间隔
        int step = IntOp::power(n, d);

        // 检查是否靠近边界
        int bound = v.nearBound(d);

        // 如果不是边界
        if (bound == 0)
        {
            // 相邻点作为 ghost cell，此时 ghost cell 真实存在
            ghost[0] = v.value(-2 * step);
            ghost[1] = v.value(-step);
            ghost[3] = v.value(step);
            ghost[4] = v.value(2 * step);
        }
        // 如果接近边界
        else
        {
            // 每一维上有两个边界条件 -1 表示靠近下边界，1 表示靠近上边界
            int sign = (bound > 0) ? 1 : -1;

            // 填充不在边界的 ghost cell
            ghost[2 - sign * 2] = v.value(-2 * sign * step);
            ghost[2 - sign] = v.value(-sign * step);

            // 到边界的距离
            int dist = 2 - IntOp::abs(bound);

            // 获得边界类型
            int type = type_[2 * d + ((bound > 0) ? 1 : 0)];

            // 距离边界为 1 的情况
            if (dist == 1)
            {
                ghost[2 + sign] = v.value(sign * step);
                ghost[2 + sign * 2] = oneGhost(v, type, d, sign, dist);
            }
            // 处于边界的情况
            else
            {
                ghost[2 + sign] = oneGhost(v, type, d, sign, dist);
                ghost[2 + sign * 2] = twoGhost(v, type, d, sign, dist);
            }
        }
    }

    /**
     * @file GhostFiller.h
     * @author xingyifan
     * @date 2022-08-30 13:23
     *
     * @description: 用于补足边界条件：<phi>i+ed = ... + 4<g>i+0.5ed
     * 非周期边界需要补充 <g>
     */

    Vector bdCond(Tensor<Dim> &v, Real t)
    {
        // 创建容器，初始化
        Vector bdvalue(5 * Dim, 0);

        // 循环每一个方向，填充 ghost cell
        for (int d = 0; d < Dim; d++)
        {
            // 获得指定方向上的边界条件
            bdCond(v, bdvalue.ptr() + 5 * d, d, t);
        }

        return bdvalue;
    }

    /**
     * @file GhostFiller.h
     * @author xingyifan
     * @date 2022-08-31 13:12
     *
     * @description: 计算指定维上的边界条件
     */

    void bdCond(Tensor<Dim> &v, Real *bdvalue, int d, Real t)
    {
        // 网格规模
        int n = v.width();

        // 分割宽度
        Real h = 1.0 / n;

        // 检查是否靠近边界
        int bound = v.nearBound(d);

        // 如果接近边界
        if (bound != 0)
        {
            // 每一维上有两个边界条件 -1 表示靠近下边界，1 表示靠近上边界
            int sign = (bound > 0) ? 1 : -1;

            // 包装边界函数
            FuncX g = [=](Vector x)
            {
                return g_(x, t);
            };

            // 计算边界条件
            Integrator<3, Dim> I;

            // 创建单位区域
            Domain<Dim> domain;

            // 获得边界类型
            int type = type_[2 * d + ((bound > 0) ? 1 : 0)];

            // 到边界的距离
            int dist = 2 - IntOp::abs(bound);

            // 记录边界值
            Real bd = 0;

            // 狄利克雷边界
            if (type == 0)
            {
                // 注意如果距离边界为 1，需要移动一个单元计算面积分
                bd = (Real)I.FaceAverage(domain, v, d, (1 + dist) * sign, g);
            }
            // 诺依曼边界
            if (type == 1)
            {
                // 注意如果距离边界为 1，需要移动一个单元计算面积分
                Real phii = (Real)I.CellAverage(domain, v, g, d, dist * sign);
                Real phii1 = (Real)I.CellAverage(domain, v, g, d, (dist - 1) * sign);
                Real phii2 = (Real)I.CellAverage(domain, v, g, d, (dist - 2) * sign);
                Real phii3 = (Real)I.CellAverage(domain, v, g, d, (dist - 3) * sign);
                Real phiface = (Real)I.FaceAverage(domain, v, d, (1 + dist) * sign, g);

                bd = (-415 * phii + 161 * phii1 - 55 * phii2 + 9 * phii3 + 300 * phiface) / (72 * h);
            }

            // 距离边界为 1 的情况
            if (dist == 1)
            {
                if (type == 0)
                {
                    bdvalue[2 + sign * 2] = bd * 4;
                }
                if (type == 1)
                {
                    bdvalue[2 + sign * 2] = bd * 6 * h / 5;
                }
            }
            // 处于边界的情况
            else
            {
                if (type == 0)
                {
                    bdvalue[2 + sign] = bd * 4;
                    bdvalue[2 + sign * 2] = bd * 16;
                }
                if (type == 1)
                {
                    bdvalue[2 + sign] = bd * 6 * h / 5;
                    bdvalue[2 + sign * 2] = bd * 6 * h;
                }
            }
        }
    }

private:
    /**
     * @file GhostFiller.h
     * @author xingyifan
     * @date 2022-08-31 22:41
     *
     * @description: oneGhost 获取 i+ed
     *
     * type 类型
     * d 方向
     * sign 符号
     * dist 到边界的距离
     */

    Real oneGhost(Tensor<Dim> &v, int type, int d, int sign, int dist)
    {
        // 存放结果
        Real res = 0;

        // 网格宽度
        int n = v.width();

        // 在此方向上移动一个点的位置间隔
        int step = IntOp::power(n, d);

        // 初始指标修正
        int i = dist * sign * step;

        // 狄利克雷条件
        if (type == 0)
        {
            res = (-13 * v.value(i) + 5 * v.value(i - sign * step) - v.value(i - 2 * sign * step)) / 3;
        }
        // 诺依曼条件
        if (type == 1)
        {
            res = (5 * v.value(i) + 9 * v.value(i - sign * step) -
                   5 * v.value(i - 2 * sign * step) + v.value(i - 3 * sign * step)) /
                  10;
        }
        // 周期边界
        if (type == 2)
        {
            // 临时存放坐指标
            int tmp = v.index(d);

            // 如果靠近上边界，取下边界单元；否则取上边界单元
            v.index(d) = (sign == 1) ? 0 : n - 1;

            res = v.value();

            // 恢复坐标
            v.index(d) = tmp;
        }

        return res;
    }

    /**
     * @file GhostFiller.h
     * @author xingyifan
     * @date 2022-08-31 22:41
     *
     * @description: twoGhost 获取 i+2ed
     */

    Real twoGhost(Tensor<Dim> &v, int type, int d, int sign, int dist)
    {
        // 存放结果
        Real res = 0;

        // 网格宽度
        int n = v.width();

        // 在此方向上移动一个点的位置间隔
        int step = IntOp::power(n, d);

        // 初始指标修正
        int i = dist * sign * step;

        // 狄利克雷条件
        if (type == 0)
        {
            res = (-70 * v.value(i) + 32 * v.value(i - sign * step) - 7 * v.value(i - 2 * sign * step)) / 3;
        }
        // 诺依曼条件
        if (type == 1)
        {
            res = (-75 * v.value(i) + 145 * v.value(i - sign * step) -
                   75 * v.value(i - 2 * sign * step) + 15 * v.value(i - 3 * sign * step)) /
                  10;
        }
        // 周期边界
        if (type == 2)
        {
            // 临时存放坐指标
            int tmp = v.index(d);

            // 如果靠近上边界，取下边界单元；否则取上边界单元
            v.index(d) = (sign == 1) ? 1 : n - 2;

            res = v.value();

            // 恢复坐标
            v.index(d) = tmp;
        }

        return res;
    }

public:
    int type_[2 * Dim]; // 边界类型 0 狄利克雷 1 诺依曼 2 周期
    FuncXt g_;          // 边界函数
};